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Soft gamma repeaters (SGRs) have unique properties that make them intriguing targets for gravi- 
tational wave (GW) searches. They are nearby, their burst emission mechanism may involve neutron 
star crust fractures and excitation of quasi-normal modes, and they burst repeatedly and sometimes 
spectacularly. A recent LIGO search for transient GW from these sources placed upper limits on a 
set of almost 200 individual SGR bursts. These limits were within the theoretically predicted range 
of some models. We present a new search strategy which builds upon the method used there by 
"stacking" potential GW signals from multiple SGR bursts. We assume that variation in the time dif- 
ference between burst electromagnetic emission and burst GW emission is small relative to the GW 
signal duration, and we time-align GW excess power time-frequency tilings containing individual 
burst triggers to their corresponding electromagnetic emissions. Using Monte Carlo simulations, we 
confirm that gains in GW energy sensitivity of N 1 ^ 2 are possible, where N is the number of stacked 
SGR bursts. Estimated sensitivities for a mock search for gravitational waves from the 2006 March 
29 storm from SGR 1900+14 are also presented, for two GW emission models, "fluence-weighted" 
and "flat" (unweighted). 

PACS numbers: 04.80.Nn, 07.05.Kf 95.85.Sz 



I. INTRODUCTION 

Soft gamma repeaters (SGRs) are promising potential 
sources of gravitational waves (GWs). They sporadi- 
cally emit brief (« 0.1s) intense bursts of soft gamma- 
rays with peak luminosities commonly up to 10 42 erg/s. 
Three of the five known galactic SGRs have produced 
rare "giant flare" events with initial bright, short (w 0.2 s) 
pulses followed by tails lasting minutes, with peak lu- 
minosities between 10 44 and 10 47 erg/s. According to 
the "magnetar" model SGRs are neutron stars with ex- 
treme magnetic fields ~ 10 15 G[1]. Bursts may result 
from the interaction of the star's magnetic field with 
its solid crust, leading to crustal deformations and oc- 
casional catastrophic cracking [2-4] with subsequent ex- 
citation of nonradial star modes [5-7] and the emission of 
GWs [6-9]. For reviews, see [10, 11]. 

The LIGO Scientific Collaboration recently completed 
a search for transient gravitational waves associated with 
almost 200 individual electromagnetic SGR triggers [12]. 
That search did not detect GW, but it explicitly targeted 
neutron star /-modes, placed the most stringent upper 
limits on transient gravitational wave amplitudes at the 
time it was published, and set isotropic emission energy 
upper limits that fell within the theoretically predicted 
range of some SGR models [7] . 

In this paper we extend that work and describe a new 
electromagnctically triggered search method for gravita- 
tional waves from multiple SGR bursts. Triggered grav- 
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itational wave searches assume that gravitational waves 
associated with an electromagnetic astrophysical event 
would reach earth at approximately the same time as 
light from the event. In addition, source sky location is 
also known. 

Knowledge of time and sky position can be a great 
advantage to gravitational wave searches [13]. It allows 
us to calculate the detector response functions, allowing 
us to estimate more relevant upper limits on GW emis- 
sion using simulated signals from the source. Also, up- 
per limits are typically lower than for untriggered all-sky 
searches, largely because they are more robust to loud 
glitches. Finally, searches can target known astrophysi- 
cal events. If the distance to the source is known, results 
can be given in terms of isotropic gravitational wave en- 
ergy emitted from the source instead of strain amplitude 
at the detector. This ties the search to the astrophysical 
source instead of the detector on Earth. All of these ad- 
vantages apply to searches for gravitational waves from 
SGR bursts. 

The new analysis method described here, "Stack-a- 
flare," builds upon the analysis pipeline ("Flare pipeline") 
used in [12] and described in [14, 15] by attempting to 
"stack" potential gravitational wave signals from multi- 
ple SGR bursts. To stack N bursts, we first generate 
N excess power time-frequency tilings. These are 2- 
dimcnsional matrices in time and frequency generated 
from the two detectors' data streams. Each tiling ele- 
ment gives an excess power estimate in the GW detector 
data stream in a small period of time St and a small 
range of frequency Sf. The time range of each tiling is 
chosen to be centered on the time of one of the target EM 
bursts in the storm. We then align these N tilings along 
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the time dimension so that times of the target EM bursts 
coincide, and perform a weighted addition according to 
a particular GW emission model. 

We hope to improve the search sensitivity by com- 
bining potential gravitational wave signals from separate 
bursts in an attempt to increase the signal-to-noise ra- 
tio, increasing the probability of detection and placing 
more stringent constraints on theoretical models via up- 
per limits. We expect that this method would be suitable 
for performing searches using data from intcrfcrometric 
detectors such as LIGO's, for gravitational waves asso- 
ciated with SGR storm events such as the 2006 March 
29 SGR 1900+14 storm. Figure 1 shows the storm light 
curve obtained by the BAT detector aboard the Swift 
satellite [16], and Figure 2 illustrates the stacking proce- 
dure with the four most energetic bursts from the storm, 
showing the main search timescales. 

This paper is organized as follows. In Section II we dis- 
cuss the general strategy of multiple SGR burst searches. 
In Section III we describe two versions of the analysis 
method ("Stack-a-flare") , both of which are built upon 
the Flare pipeline. One version coherently stacks GW 
time series associated with electromagnetic bursts in the 
storm, while the other stacks incoherently. We character- 
ize the two versions using simulations in gaussian noise, 
demonstrating the strengths and weaknesses of each and 
showing that relatively weak signals which could not be 
detected in the individual burst search can easily be de- 
tected with the new method. Gains in GW energy sen- 
sitivity of N 1 / 2 are feasible with the incoherent version 
even in the presence of realistic relative timing uncer- 
tainties between SGR bursts, where N is the number of 
stacked SGR bursts. Finally, in Section IV we examine 
the SGR 1900+14 storm light curve of Figure 1 in detail 
and present estimated search sensitivities for a simulated 
search for gravitational waves from the storm, for two 
GW emission models: "fluence-wcighted" and "flat" (un- 
weighted) . 

II. STRATEGY 

Before discussing Stack-a-flare method implementation 
details, in this section we discuss the strategy of a gen- 
eralized multiple SGR GW search. Due to similarities 
between the individual burst and multiple burst searches 
in terms of astrophysics, goals, and implementation, the 
search choices will be similar to those followed in [12]. 

A. Parameter space of the target GW signals 

As in the individual SGR search, multiple SGR 
searches could target neutron star fundamental mode 
ringdowns (RDs) predicted in [5-7, 17] as well as un- 
known short-duration GW signals. The former are inter- 
esting targets because fundamental modes are the modes 
most efficiently damped by GW emission. The latter 



are interesting targets because the flare mechanism may 
involve some other violent motion inside the star, for ex- 
ample matter carried along with magnetic field lines as 
they rearrange themselves inside the crust. We assume 
that given a neutron star, /-mode frequencies and damp- 
ing timescales would be similar from event to event, and 
that unknown signals would at least have similar central 
frequencies and durations from event to event. 

We thus focus on two distinct regions in the target 
signal time-frequency parameter space. The first region 
targets ~(100-400)ms duration signals in the (1-3) kHz 
band, which includes /-mode RD signals predicted in [18] 
for ten realistic neutron star equations of state. We 
choose a search band of (1-3) kHz for RD searches, with 
a 250 ms time window which was found to give optimal 
search sensitivity [15]. The second region targets ~(5- 
200) ms duration signals in the (100-1000) Hz band. The 
target durations are set by prompt SGR burst timescales 
(5 ms to 200 ms) and the target frequencies are set by 
the detector's sensitive region. We search in two bands: 
(100-200) Hz (probing the region in which the detectors 
arc most sensitive) and (100-1000) Hz (for full spectral 
coverage below the RD search band) using a 125 ms time 
window. In all, we could search in three frequency bands. 

We could explore these areas in the parameter space 
using the same twelve simulated waveform types used for 
setting upper limits in the individual SGR burst search, 
described in [12, 15]: linearly and circularly polarized 
RDs with t = 200 ms and frequencies in the range 1- 
3 kHz; and band- and time-limited white noise bursts 
(WNBs) with durations of 1 1 ms and 100 ms and fre- 
quency bands matched to the two low frequency search 
bands. Polarization angle will be chosen randomly for 
every compound injection in initial searches. 

It seems plausible to assume that, for a given neutron 
star, /-mode frequencies and damping timescales would 
be similar from event to event. For pure /-modes these 
quantities depend only on the global structure of the star, 
essentially the mean density [5], and the magnetic cor- 
rections depend on integrals of the field over the entire 
stellar interior (e.g. [19]). However, the major motiva- 
tion for the low frequency unknown-waveform portion 
of the parameter space is stochastic gravitational wave 
emission arising from violent events in the neutron star 
crust. Therefore, we will not assume similar waveforms 
from event to event in the unknown-waveform portion, 
although we will assume similar central frequencies and 
durations. 



B. On-source region 

In a search, we divide the GW data into an on-source 
time region, in which GWs associated with a given burst 
could be expected, and a background time region. On- 
source and background segments are analyzed identically 
resulting in lists of analysis events. 

For the individual SGR search, trigger timing precision 
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FIG. 1: BAT light curve of the SGR 1900+14 storm event, 1ms bins, 15-100 keV. The light curve shows the BAT event 
data, from sequence 00203127000, from approximately 20 s after the start of the sequence to its end. The data are publicly 
available [20] . Several of the largest bursts, durations greater than 500 ms, are considered intermediate flares. The x-axis 
times are relative to 2006-03-29 02:53:09. 9±0. 5 s UT at the satellite. Times used in the analysis itself have been corrected for 
time-of-fiight delays to the geocenter, which change by 0.35 ms over the ~30.5s course of the light curve. 
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FIG. 2: Individual EM bursts inform the stacking of GW data. This figure suggests the stacking procedure and explicitly shows 
search timescales. The top four plots are EM light curve time series around individual bursts beginning at 2.0 s, 16.6 s, 18.3 s, 
and 22.3 s in Fig. 11. Simulated GW ringdowns in the fluence-weighted model are superposed. The bottom plot shows the EM 
time series simultaneously, and the sum of the hypothetical GW signals. The on-source region of ±2s is shaded. In a search, 
GW data corresponding to the EM time series are transformed to time-frequency power tilings before being added together 
and therefore there is no dependence on phase-coherence of GW signals in the analysis; this transformation is not illustrated. 



on the order of a second was handled with 4 s on-source 
regions. For a multiple SGR search, significantly higher 
precision in relative trigger times between burst events 
in the stack will be required. However, a common bias in 
trigger times shared by all bursts in the stacking set (ab- 
solute timing) can be handled with an adequately large 
(e.g. 4 s) on-source region. 

In general, imprecision in trigger times comes princi- 
pally from two sources: satellite to geocenter light cross- 
ing delay and arbitrariness of the satellite trigger point 
in the light curve. It is possible to decrease both un- 
certainties through careful analysis of satellite data. If 
necessary, light crossing times at satellites can be prop- 
agated to the geocenter (and subsequently to any given 



interferometer) using the appropriate cphcmcris. If satel- 
lite data is public, light curve analysis can yield a specific 
timing point, for example the start of the steep burst rise. 
Increased absolute timing precision could allow us to use 
smaller on-source regions with durations set by theoret- 
ical predictions of time delay between electromagnetic 
and gravitational wave emission from SGR bursts. How- 
ever, this could exclude astrophysical models with larger 
timing delays, for instance if the initial excitation trans- 
fers energy to the /-mode through some slow and inef- 
ficient mode-mode coupling due to selection rules. And 
it turns out there is little benefit to be gained: We have 
performed Monte Carlo simulations comparing ±2 s and 
±ls on-source regions. Reducing the on-source region 
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from 4 s to 2 s resulted in only a 2% reduction in ampli- 
tude upper limits, on average over 24 trials with various 
waveform types. 

We make one new assumption about the nature of indi- 
vidual bursts of gravitational waves from SGRs: that the 
time delay between each gamma-ray burst and the associ- 
ated GW burst does not vary by more than about 100 ms, 
less than duration of the shortest coherently modeled sig- 
nal and comparable to the star's Alfven crossing time. A 
variation this large would severely degrade detection ef- 
ficiency in the case of 11ms duration WNBs, but as we 
shall see in Section III F, detection of neutron star RDs 
(with r <~ 200 s) remains robust. 

C. Background region 

As with the individual burst search, the background 
region serves three purposes [14, 15]: 

1. it is used to estimate statistics of the power tiling as 
a function of frequency for use in the Flare pipeline; 

2. it provides false alarm rate (FAR) estimates from 
which the significance of the loudest on-source anal- 
ysis event can be determined; 

3. it is representative noise into which simulated wave- 
forms can be injected for estimating upper limits. 

In the course of validating the individual SGR search 
we showed that 1000 s of data on either side of an on- 
source region produce sufficient estimates of the power 
tiling statistics [15]. This requirement and the estimation 
procedure are unchanged in the multiple SGR search, so 
±1000 s of data will again suffice for this purpose. The 
background region required for injecting simulations to 
estimate upper limits may depend on the system be- 
ing modeled and the desired statistical precision; for the 
mock SGR 1900+14 storm search we describe below, 
±1000 s of background is sufficient. The background re- 
gion required for FAR estimates depends on the desired 
precision of FAR estimates. Estimating the FAR of a 
very large on-source analysis event requires a larger back- 
ground than estimating the FAR of a small on-source 
analysis event, for a given level of precision. 

D. GW emission models 

In a multiple burst search we must choose which bursts 
to include in the set and how to weight them. Since GW 
energy £qw is unknown we are left to attempt to choose 
weightings via other observables which we hope will cor- 
relate to it. As with the individual burst search, we as- 
sume that the SGR burst sample is comprised of bursts 
occurring within some specified time range defined by the 
observatory's science run schedule, and not necessarily 
from the same SGR source. We will refer to a set of SGR 



bursts to be included in the multiple burst search, along 
with a weighting strategy, as a GW emission model. 

We could use Occam's razor to select GW emission 
models such as the following: 

si. flat model — use every detected and confirmed burst 
from a given SGR source within the time range, with 
equal weighting. In practice we need to choose a cut- 
off on which bursts to include; 

s2. inclusive model — use every detected and confirmed 
burst within the time range, from any SGR source, 
with equal weighting; 

s3. fluence-weightcd model — use every detected and 
confirmed burst from a given SGR source within the 
time range, weighted proportional to fluence; 

s4. use a subset of component bursts from a multi- 
episodic burst event such as the SGR 1900+14 storm, 
with some weighting scheme. 

We note that the inclusive model s2 could benefit from 
a search method that was insensitive to variations be- 
tween SGR sources. For example, we would expect two 
different sources to have /-modes at different frequen- 
cies, which may not brighten corresponding pixels in a 
time-frequency tiling. We do not attempt to solve this 
problem here. 

We can also imagine stacking GW emission models 
based on specific theoretical models. A theory may pre- 
dict that there is no correlation between Eem and Eqw- 
Such a prediction corresponds to the flat model and is 
tested below. Or a theory may predict that there is 
one emission mechanism with a constant mechanical ef- 
ficiency, and thus Eem is always proportional to Eqw- 
This is treated in the fluence-weighted model and tested 
below. A theory may predict that £qw is some more 
complicated function of -Eem, the simplest example being 
a step function — for instance, a theory that small flares 
originate from magnetic reconnections far out in the mag- 
netosphere and only large flares from reconnections in or 
near the star itself [21]. However the Swift/BAT spectra 
of the SGR 1900+14 storm indicate that the neutron star 
surface participates in all flares [22] , and therefore we do 
not address such a model below (although it might be 
interesting at a later date). A theory may also predict 
that the time delay between electromagnetic and gravi- 
tational emission varies from burst to burst. Some varia- 
tion, up to the target signal durations of tens or hundreds 
of milliseconds, is accounted for in the searches described 
below. Larger variations could potentially be treated by 
sweeping over some range of time delays for each burst, 
a more complicated search which we do not yet address. 

We will neglect for the time being complicating factors 
such as: multiple injections of energy into a single burst, 
with possible correlation in gravitational wave emission 
energy; qualitatively different gravitational wave emis- 
sion in the case of intermediate flares and common bursts; 
and beaming issues. We also neglect other possibilities 
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such as correlating Eqw with flare peak luminosities or 
rise times. 



III. ANALYSIS METHOD 

Both versions of the Stack-a-flare pipeline, "T-Stack" 
and "P-Stack," consist of an extension to the Flare 
pipeline [14, 15], a simple but effective search pipeline 
based on the excess power detection statistic of [23]. 



A. Flare pipeline 

The Flare pipeline can process streams of data from ei- 
ther one or two GW detectors. First data is conditioned, 
then excess power tilings are created, and finally analy- 
sis events are constructed via clustering elements in the 
tilings. 

Data conditioning consists of zero-phase digital filter- 
ing in the time domain [24], first with a bandpass filter 
and then with a composite notch filter. The raw cali- 
brated LIGO power spectrum is colored, and is charac- 
terized by a sensitive region between ~60 Hz to ^2 kHz 
which includes a forest of narrow lines, with increasingly 
loud noise on cither side of the sensitive band. Search 
sensitivity is increased by removing these insensitive re- 
gions from the data, which would otherwise dominate 
weak signals and destroy bandwidth after transformation 
to frequency domain. We remove narrow lines associated 
with the power line harmonics at multiples of 60 Hz, the 
violin modes of the mirror suspension wires, calibration 
lines, and persistent narrow band noise sources of un- 
known origin. 

Time-frequency spectrograms arc then created from 
conditioned data for individual detectors from a series of 
Blackman- windowed discrete Fourier transforms, of time 
length St set by the target signal duration. A tile is an 
estimate of the short-time Fourier transform of the data 
at a specific time and frequency. Each column in the 
tiling corresponds to a time bin of width 5t and each 
row corresponds to a frequency bin of width Sf, both 
linearly spaced, with Sfdt — 1. Adjacent time bins over- 
lap by 0.9St to guard against mismatch between prospec- 
tive signals and tiling time bins. Larger overlaps require 
more computation and do not noticeably improve sensi- 
tivity [15]. 

In a one-detector search, we then have a complex- 
valued time- frequency tiling from which we calculate the 
real-valued one-sided PSD for every time bin. To do this 
we multiply each tile value by its complex conjugate and 
normalize the result to account for sampling frequency 
and windowing function. We discard frequency bins out- 
side of the chosen search band. 

In a two-detector search, we have two complex time- 
frequency tilings (one for each detector) from which we 
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where T represents a tiling matrix and t and / arc time 
and frequency bin indices, and (1) and (2) denote the 
detector. Here At is the gravitational wave crossing time 
difference between detectors; this term takes care of ap- 
plying the appropriate time difference between detector 
data streams in the Fourier space, with the advantage 
of permitting sub-sample time delays, which significantly 
increases the sensitivity at higher frequencies. The real 
part is kept, and normalization is applied as in the one- 
detector case. To obtain a positive-definite statistic we 
take the absolute value of each tile; this allows sensitivity 
to both strongly correlated and strongly anti-correlated 
signals in the two (potentially misaligned) detectors. 

Next, we use off-source data to remove the background 
noise power from each clement of the PSD time-frequency 
tiling. The elements are fit to a gamma distribution [15], 
and outliers above a threshold (typically four standard 
deviations) are discarded. The resulting estimate on the 
mean is subtracted from each element of the correspond- 
ing frequency bin in the PSD matrix, giving a matrix of 
excess power (or "cross excess power" in the two-detector 
case) . 

To create analysis events we use a simple 2-dimensional 
clustering algorithm, allowing retention of signal energy 
which might otherwise be fragmented in the case of ex- 
tended signals in the time- frequency plane. Analysis 
events correspond to discrete clusters found by the al- 
gorithm, and include information on cluster central fre- 
quency, central time, bandwidth, duration, and so forth. 
The detection statistic (event loudness) is the sum over 
the cluster of tile significance. 



B. T-Stack-a-fiare pipeline 

The T-Stack pipeline combines burst events in the time 
domain. 

For each of N burst events a trigger time is determined. 
For a given GW detector, ./V time scries segments con- 
taining those trigger times are then aligned to the trigger 
times, weighted according to antenna factor, and added 
together. The resulting summed time series (either one 
or two, depending on how many detectors are included 
in the search) are then fed to the Flare pipeline. 

As will be described below, the T-Stack pipeline has 
the advantage of achieving higher sensitivity in gaussian 
noise, but the disadvantage of being extremely sensitive 
to timing inaccuracies and variations in GW signal wave- 
form from burst to burst. 

This makes it a potentially viable choice for analyz- 
ing multi-episodic events — in which a single contiguous 
100 /xs-binned light curve might provide adequate timing 
precision — but a poor choice for analyzing isolated burst 
events or incoherent signals such as band-limited WNBs. 
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FIG. 3: Diagram of the T-Stack version of the Stack-a-fiare pipeline. The T-Stack pipeline has a thin layer added before the 
Flare pipeline in which gravitational wave data time series containing SGR burst event triggers are aligned on the trigger times 
and added together. These stacked time series are made for each detector and then run through the Flare pipeline as normal. 



C. P-Stack-a-flare pipeline 

The P-Stack pipeline combines burst events in the fre- 
quency domain. 

We determine a trigger time for each of N burst events 
based on gamma-ray data. Each of N time series con- 
taining those triggers is processed with the Flare pipeline, 
up to the clustering algorithm, exactly as in an individ- 
ual SGR burst search. Antenna factors are applied at this 
time. The result is N time-frequency significance tilings, 
with different noise but hopefully with similar or identi- 
cal signals. The N significance tilings are then aligned 
to the trigger time and added together. The combined 
significance tiling is then fed through the Flare pipeline 
clustering algorithm with a fixed fraction of tiles to in- 
clude in the clustering (e.g. 0.1%). A fixed fraction of 
tiles is used instead of a fixed loudness threshold value 
because the variance of the tile loudness distribution at 
a given frequency increases with N. 

As will be described below, the P-Stack pipeline has 
the advantage of being relatively insensitive to timing in- 
accuracies or differences in waveform from burst to burst, 
but it has less sensitivity than the T-Stack pipeline for the 
(possibly unrealistic) precisely-known timing case, with 
deterministic waveforms. 



D. Loudest event upper limits 

As in the individual SGR search, in the absence of 
a detection we can estimate loudest event upper limits 
[25] on GW root-sum-squared strain h rss incident at the 
detector, and GW energy emitted isotropically from the 
source assuming a nominal source distance. 

We estimate loudest-event upper limits [25] on GW 
root-sum-squared strain h TSS incident at the detector. We 
can construct simulations of impinging GW with a given 
h Iss . Following [26] 

^rss = ^rss+ + ^rssx ' (2) 



where e.g. 

/oo 
\h + \ 2 dt (3) 
-oo 

and /i+ jX (i) are the two GW polarizations. The rela- 
tionship between the GW polarizations and the detector 
response h(t) to an impinging GW from a polar angle 
and azimuth ((?, cj>) and with polarization angle ip is: 

h{t) = F + {6, 0, ip)h+(t) + F x (6, <f>, i>)h x (t) (4) 

where F+(0,<j),ip) and F x {9 1 4>,ip) are the antenna func- 
tions for the source at (8, 0) [27]. At the time of the 
storm, the polarization-independent RMS antenna re- 
sponse (i* 1 ? + F^) 1 / 2 , which indicates the average sensi- 
tivity to a given sky location, was 0.39 for LIGO Hanford 
observatory and 0.46 for the LIGO Livingston observa- 
tory. 

We can also set upper limits on the emitted isotropic 
GW emission energy £"gw a t a source distance R associ- 
ated with h+(t) and h x (t) via [28] 

r 3 r°° / . \ 
E GW = 4TrR 2 —J ((h+) 2 + (h x ) 2 ) dt. (5) 

The upper limit is computed in a frcqucntist frame- 
work following the commonly used procedure of injecting 
simulated signals in the background data and recovering 
them using the search pipeline (see for example [29, 30]). 

An analysis event is associated with each injected sim- 
ulation, and compared to the loudest on-source analysis 
event. The GW strain or isotropic energy at e.g. 90% 
detection efficiency is the strain or isotropic energy at 
which 90% of injected simulations have associated events 
louder than the loudest on-source event. 

"Efficiency curves" are constructed by the Flare 
pipeline by repeatedly analyzing 4 s background seg- 
ments, each containing a single simulation created with a 
range of values, and comparing the loudest simula- 
tion analysis event within 100 ms (for RDs) or 50 ms (for 
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FIG. 4: Diagram of the P-Stack version of the Stack-a-fiare pipeline. The P-Stack pipeline has a thin layer added after the Flare 
pipeline in which gravitational wave data significance tilings containing SGR burst event triggers are aligned on the trigger 
times and added together. Stacked significance tilings can then be run through the Flare pipeline clustering algorithm. 
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FIG. 5: Example efficiency curve generated for the Monte 
Carlo experiment investigating Stack-a-flare sensitivity vs. TV. 
This example curve was made with the T-Stack pipeline using 
compound injections of N — 200 f 590 Hz circularly polarized 
ringdowns injected into gaussian noise with a = 1. Each 
efficiency curve was constructed using 20 amplitude scaling 
factors and 20 trials at each h las amplitude. Binomial error 
bars are given by \fr{l — r)/M where r is the fraction de- 
tected at a given injection amplitude and M is the number of 
injections in the Monte Carlo. 



WNBs) of the known injection time to the loudest on- 
source analysis event [15]. An example efficiency curve is 
shown in Figure 5. The range of h ISS values must be cho- 
sen so that the smallest value produces simulations that 



are always lost in the noise, and the largest value pro- 
duces simulations that are typically detected with large 
signal-to-noise ratios. The Eq^ or /r,^™value at 90% de- 
tection efficiency (E™$ or h^°) occurs where 90% of 
the loudest simulation analysis events are larger than the 
loudest on-source event. 

For any given on-source region this results in four ar- 
rays of numbers, each of which has length equal to the 
number of injected simulations used to estimate the up- 
per limit. The first contains the h ISS values of injected 
simulations. The second contains the calculated £gw 
values of injected simulations, or Eq^/R 2 if the distance 
to the source R is not known. The third contains the 
loudness of the analysis event associated with the injected 
simulation. The fourth contains boolean values indicat- 
ing whether the associated analysis event was larger then 
the loudest on-source analysis event or not. 

The hrss and boolean (or the -Egw & n d boolean) ar- 
rays can be used to construct the efficiency curve, with 
the h ISS (or -Eqw) values on the x-axis. The y-axis in- 
dicates the fraction of analysis events associated with an 
injected simulation of h rss as given by the x-axis which 
are louder than the loudest on-source event. In the case 
of simulation h vss values which range over a discrete set 
of scale factors, the y-axis value is simply this fraction. 
Binomial error bars may be added to these data points. 

However, we are typically interested in the or 
^rss % value, that is the h rss of simulations whose associ- 
ated analysis event is louder than the loudest on-source 
analysis event 90% or 50% of the time. Since we don't 
know this value ahead of time, it is necessary to inter- 
polate between the h ras values associated with the dis- 
crete scale factors; we do this by fitting with a sigmoid 
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function. The Flare efficiency curve fitting routine uses 
two functions to perform these fits: a four-parameter fit 
based on the logistics function, and a five-parameter fit 
based on the complementary error function. The models 
were chosen on empirical grounds, and details are given 
in [15]. 

We can follow the same procedure for the multiple 
burst search. The only difference is the need to mea- 
sure the h TSS or -Eqw of a compound injection, instead of 
a simple (single) injection. 



E. Sensitivity dependence on N 

The matched filter amplitude signal-to-noise ratio 
(SNR) is defined in the frequency domain as [31] 



Kf? 

Sn(f) 



df 



1/2 



(6) 



where h(f) is the Fourier transform of the signal time se- 
ries and S n (f) is the noise power spectral density. Here, 
the numerator is the square root of the power in the sig- 
nal. In gaussian noise with zero mean, S n (f) — a 2 , a con- 
stant. Since the standard deviation a of gaussian noise 
goes as the square root of N and the amplitude of identi- 
cal stacked signals goes as N, we expect the SNR of the 
optimal T-Stack algorithm for the recovery of identical 
signals from noise to go as TV 1 / 2 . 

Whereas the T-Stack pipeline stacks amplitude, the 
P-Stack pipeline stacks power. The background tiles in 
the power tiling at each individual frequency bin can be 
modeled as Gamma-distributed noise [15], for which the 
variance also goes as N, so we expect the power signal-to- 
noise ratio to increase as N 1 / 2 . Since the amplitude goes 
as the square root of the power, we expect the P-Stack 
amplitude sensitivity to increase as jV 1 / 4 . 

We tested these predictions by injecting N stacked 
1590 Hz 200 ms r ringdowns into gaussian noise with 
a = 1. We then constructed efficiency curves determin- 
ing the injection h TSS at 50% and 90% detection efficiency. 
Each efficiency curve was constructed using 20 amplitude 
scaling factors and 20 trials at each h rss amplitude. An 
example efficiency curve is shown in Figure 5. We then 
fit the 50% and 90% detection efficiency level results as 
functions of iV to a two-parameter power law of the form 
y = AN B . The results for both the T-Stack and P- 
Stack pipelines are shown in Figure 6 at 90% detection 
efficiency. The fit for the T-Stack pipeline gives a sen- 
sitivity dependence in amplitude at both detection effi- 
ciency levels of nearly N 1 / 2 (N 0A9 and TV 55 for 50% 
and 90% detection efficiencies respectively), confirming 
our prediction. This corresponds to an improvement in 
energy of a factor of N. The fit for the P-Stack pipeline 
gives a sensitivity dependence in amplitude at both de- 
tection efficiency levels of nearly TV 1 / 4 (TV 24 and N - 27 
for 50% and 90% detection efficiencies respectively) , con- 
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FIG. 6: T-Stack and P-Stack sensitivity dependence on N 
at 90% detection efficiency, for 1590 Hz r=200ms ringdowns 
in gaussian noise with a = 1. The T-Stack pipeline show 
a sensitivity dependence of nearly N 1 ^ 2 and the results for 
the P-Stack pipeline show a sensitivity dependence of nearly 
N 1/4 . All fits excluded the point N = 1. Results at 50% 
detection efficiency look similar. 



firming our prediction. This corresponds to an improve- 
ment in energy of a factor of N 1 / 2 . 

We repeated the experiment for 100 ms duration 100- 
1000 Hz band-limited WNBs. In this case, we expected 
the coherent T-Stack pipeline to undcrperform the the P- 
Stack pipeline on these independently-generated stochas- 
tic incoherent signals. As expected, we found that the 
T-Stack pipeline shows no improvement as N increases, 
while the P-Stack pipeline show the same TV 1 / 4 sensitiv- 
ity dependence seen in the coherent ringdown case (TV ' 36 
and N - 23 for 50% and 90% detection efficiencies respec- 
tively) . The results are shown in Figure 7; they illustrate 
the relative model-independence of the P-Stack pipeline. 



F. Sensitivity dependence on timing errors 

The T-Stack pipeline attains optimal sensitivity gains 
with increasing N because it performs a phase coherent 
addition of signals. We have shown that the P-Stack 
pipeline attains its N 1 / 2 energy sensitivity performance 
even in the case of stacked signals that are not coherent 
such as independently-generated white noise bursts. For 
identical signals such as ringdowns, errors in the relative 
times between stacked signals cause breakdown of phase 
coherence. 

We performed Monte Carlo simulations using a sim- 
ulated burst series with N — 20 equal-amplitude ring- 
downs, and allowing the timing error to range. We char- 
acterized 1090 Hz t = 200 ms and 2590 Hz r = 200 ms 
circularly polarized ringdowns, corresponding to the low 
and high frequency ranges in the signal parameter space. 
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FIG. 7: T-Stack and P-Stack sensitivity dependence on TV, 
at 90% detection efficiency, for 100-1000 Hz 100 ms duration 
white noise bursts in gaussian noise with a = 1. The re- 
sults for the T-Stack pipeline show a sensitivity dependence 
of nearly TV (flat dependence), and the results for the P-Stack 
pipeline show a sensitivity dependence of nearly TV 1 / 4 , as in 
the coherent ringdown case. All fits exclude the point N — 1. 
Results at 50% detection efficiency look similar. 



Timing errors were randomly chosen for each ringdown 
from a normal distribution with a ranging from 10 fj,s to 
100 ms, and were applied as a timing shift to the given 
ringdown. At each timing uncertainty we constructed ef- 
ficiency curves using 20 amplitude scaling factors and 20 
trials at each h ISS amplitude. 

The results are displayed in Figure 8 (1090 Hz ring- 
downs) and Figure 9 (2590 Hz ringdowns) at 90% detec- 
tion efficiency. As expected, the P-Stack method is inde- 
pendent of timing error, up until large timing errors on 
the order of the signal duration. The T-Stack pipeline, 
on the other hand, shows a pronounced dependence on 
timing error, especially in the case of high frequency sim- 
ulations. Each plot shows data for both T-Stack and 
P-Stack pipelines, and finds the equal-sensitivity timing 
error (P-Stack and T-Stack curve intersection point) us- 
ing polynomial fits. 

For the T-Stack pipeline to be effective at 1090 Hz, 
apparently, timing error must be < 100 fis at one- sigma. 
For the T-Stack pipeline to be effective at 2590 Hz, timing 
error must be < 50 (is at one-sigma. For the TV = 20 case 
shown, with no timing error the T-Stack pipeline is about 
1.5 times more sensitive than the P-Stack pipeline. 



G. Implications of the pipeline characterization 

We summarize the implications from characterizing the 
T-Stack and P-Stack pipelines. We envision four possible 
types of stacked SGR searches: 

1. High frequency (1000-3000 Hz) searches for ringdown 
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FIG. 8: T-Stack and P-Stack sensitivity versus timing error, 
for 1090 Hz t = 200 ms circularly polarized RD, TV = 20, 



for h 



90% 



T-Stack is more sensitive for small timing errors, 



but degrades. The crossover point is noted; for T-Stack to 
be effective at 1090 Hz, apparently timing error must be < 
100 /is at one-sigma. T-Stack results level off at high timing 
errors (greater than ~ 2 x 10~ 4 , or ~90 degrees of phase) 
because the Monte Carlo effectively randomizes the phases of 
the stacked signals. The green region shows the approximate 
range of timing uncertainties, shown in Table I. Results at 
50% detection efficiency look similar. 



burst emission, for single SGR storm events (ringdown 
upper limits); 

2. Low frequency (100-1000 Hz) searches for stochastic 
burst emission, for single SGR storm events (band- 
and time-limited WNB upper limits); 

3. High frequency (1000-3000 Hz) searches for ringdown 
burst emission, for isolated, time-separated SGR 
bursts (ringdown upper limits); 

4. Low frequency (100-1000 Hz) searches for stochas- 
tic burst emission, for isolated, time-separated SGR 
bursts (band- and time- limited WNB upper limits). 

We distinguish storm events from isolated events be- 
cause it is much easier to get precise relative times for 
storm events. 

We have found that the P-Stack pipeline should be ef- 
fective in any of these cases, with an energy sensitivity 
gain over the individual burst search of approximately 
TV 1 / 2 . The T-Stack pipeline shows an energy sensitivity 
improvement of approximately TV, but only if the target 
signals are coherent, and only if the relative timing be- 
tween SGR GW burst events can be known to high preci- 
sion (<100^s at 1090Hz and <50 /is at 2590Hz). As we 
will sec, this timing precision requirement is stringent. In 
the future we might search for a method of overcoming 
the T-Stack timing precision limitation, possibly by us- 
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FIG. 9: T-Stack and P-Stack sensitivity versus timing error, 
for 2590 Hz r = 200 ms circularly polarized RD, TV = 20, 



for h 



T-Stack is more sensitive for small timing errors, 



but degrades. The crossover point is noted; for T-Stack to 
be effective at 2590 Hz, apparently, timing error must be < 
50 /xs at one-sigma. T-Stack results level off at high timing 
errors (greater than ~ 1 x 10~ 4 , or ~90 degrees of phase) 
because the Monte Carlo effectively randomizes the phases of 
the stacked signals. The green region shows the approximate 
range of timing uncertainties, shown in Table I. Results at 
50% detection efficiency look similar. 



ing additional computational resources to perform time 
shifts between the time series in the stack. 



TABLE I: The 18 most energetic peaks in the storm light 
curve included in the analysis, ordered by time. First column 
is peak number. Second column is the time in seconds rela- 
tive to the beginning of the light curve shown (see Figure 1), 
with the one-sigma uncertainty from the rising edge fit given. 
Third column is light travel time from satellite to geocenter in 
ms, which is applied before analysis. Fourth column is inte- 
grated burst counts. Fifth column is estimated burst fluence, 
based on a conservative estimate of 1.0 x 10 _4 ergcm -2 for 
the total fluence of the storm by the Konus-Wind team, in the 
20-200 keV range [32] . Sixth column is fitted rising edge slope 
in counts/s. Seventh column is the approximate duration of 
the peak in seconds. 
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IV. SGR 1900+14 STORM MOCK SEARCH 

In this section we describe a mock stacking search with 
the P-Stack pipeline, using simulated LIGO data, for GW 
associated with the 2006 March 29 SGR 1900+14 storm 
event. We first describe a procedure for estimating burst 
start times and fluences from the light curve. We then 
discuss GW emission models. We finally present sensi- 
tivity estimates. 

A. Light curve and GW emission models 

Data from the BAT detector on the Swift satellite are 
publicly available. In Figure 1, we show the storm light 
curve in photon counts in the 15-100 keV band with 1 ms 
time bins. SGR bursts lasting longer than 500 ms are 
typically considered "intermediate flares"; the storm con- 
tained both intermediate flares and common bursts. 

Before choosing GW emission models for stacking, we 
gathered information from the light curve, most impor- 
tantly consistent relative burst times. Though the BAT 
timing resolution is 100 fis, additional resolution would 
not improve our estimates of relative burst times as we 



shall see. We also measured integrated counts under each 
burst, a quantity roughly proportional to fluence. These 
efforts were complicated by overlapping and non-uniform 
bursts in the noisy light curve. We found that a 30- 
samplc running average of the 1 ms-binned light curve 
aided aspects of the analysis. 

We chose the point of intersection of the rapid rising 
edges of each burst with the light curve noise floor as 
plausible and consistent estimates of the times at which 
possible progenitor catastrophic neutron star events oc- 
curred. So long as GW emission is shifted less than 1.5 s 
before or after this time, it should be well within the 
on-source region. We first used a generic peak-finding 
routine to fit and find approximate peak locations; the 
smallest peaks were ignored as insignificant given the 
GW emission stacking models we ultimately chose. Lo- 
cal maxima were found by exploring the neighborhoods 
around the fitted peaks. Walking left from local maxima, 
a tuned first derivative threshold was used to determine 
the top and bottom of each burst's rising edge. The ris- 
ing edge between these bounds was then fit with a line. 
Figure 1 1 shows the fits on all bursts considered for inclu- 
sion in the analysis. We then corrected each burst time 
for the light travel time to the geocenter using the known 
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SGR 1900+14 sky position and known satellite positions. 
The light travel delay at the beginning of the light curve 
is 17.48 ms (arriving at the satellite first) and at the end 
is 17.12 ms, a change of 0.35 ms over a course of 30.5 s. A 
histogram of one-sigma uncertainty estimates is shown in 
Figure 12. These uncertainty estimates were applied in- 
dividually to respective simulations in the Monte Carlo. 

Integrated counts under each burst were also esti- 
mated, by walking to the right along the smoothed light 
curve from peak markers until nearly reaching the noise 
floor or encountering the beginning of the next burst's ris- 
ing edge. The light curve with boundaries used in the in- 
tegrated counts estimate indicated is shown in Figure 10. 

Burst integrated count uncertainties were conserva- 
tively though qualitatively estimated to be less than 5%, 
except for one overlapping burst with a burst rise at 
about 4.2 s. Uncertainty in the overlapping peak is as- 
trophysical in nature and no larger than 3% of the total 
storm A§uence. The Konus-Wind team reported a con- 
servative storm total fluence of 1 — 2 x 10~ 4 ergcm~ 2 in 
the 20-200 keV range [32] . To convert integrated counts 
to fluences, we used the most conservative estimate of 
1.0— 4crgcm~ 2 and the total integrated counts in the 
entire storm to obtain fluence per count and thus flu- 
ences for each burst in the storm. In doing so we as- 
sume bursts exhibit spectral uniformity. We wish to 
know the absolute fluence in order to set upper limits on 
7 = I E EM . Uncertainty in 7 from integrated counts 
estimates is negligible compared to the uncertainty from 
the conversion, and so we simply absorb it there. Burst 
measurements are given in Table I. 

We chose to explore two GW emission stacking models 
in this mock search: an N = 11 flat model which sets GW 
energy emission constant; and a fluence- weighted model 
comprised of the 18 brightest bursts in the light curve 
which sets 7 to be constant. The N = 11 cutoff in the 
flat model (making it a step function) is motivated by 
burst integrated counts (Figure 13). The N = 18 cutoff 
in the fluence-weighted model was motivated by Figure 
14. Including the most energetic 18 bursts accounts for 
95% of the counts in the 42 bursts considered for the 
analysis. After the 15 most energetic bursts, each ad- 
ditional burst (when ordered by energy) contributes less 
than 1% to the total. In Monte Carlo simulations we 
observed only a 3% difference (averaged over amplitude 
sensitivity estimates set via the 12 injection waveforms) 
between fluence-weighted models with cutoffs at N = 17 
and N = 35. 

We chose a fluence-weighted GW emission model to ex- 
plore the constant-7 hypothesis. We chose not to pursue 
GW emission models weighted to e.g. burst peak heights 
or burst rising edge slopes. The cross-correlation between 
estimated burst peak counts and burst integrated counts 
is 0.75 in a normalization where auto-correlations are 
unity. In the fluence-weighted emission model, assuming 
that GW emission energy is proportional to fluence we 
weight compound simulations with the square root of in- 
tegrated counts, and then normalize so that the weight of 



the most energetic burst is unity. We weight significance 
tilings on the other hand according to burst integrated 
counts before stacking them. 



B. Results 

In this section we present the results of a mock search 
with the P-Stack pipeline for GW associated with the 
2006 March 29 SGR 1900+14 storm, using simulated 
data. At the time of the storm, all three LIGO detec- 
tors were taking science-quality data. Simulated data 
modeled on data from the two LIGO 4 km detectors were 
created from gaussian noise by matching power spectra 
with the LIGO data in the frequency domain. Therefore, 
the sensitivity estimates should be similar to upper limit 
estimates from real data. 

In Figure 15 we show example cumulative histograms 
showing false alarm rates versus analysis event loudness 
for the background and the stacked on-source region. 
There are three such plots for each emission model, one 
per search band. Since the stacked livctime is 4 s here, the 
loudest on-source event occurs once per 4s, and is plot- 
ted at a y-value of 0.25 Hz. We can estimate the FAR of 
this loudest on-source analysis from the background. 

If only one emission model yields a significant event, 
or if the case is marginal, we can gain additional infor- 
mation by making a joint statement of probability. We 
can examine the background probability density function 
in a 2-dimensional space comprised of duples (A, B) of 
loudest events from corresponding background segments 
analyzed under emission models A and B, and deter- 
mine constant probability contours c = a A + B, where 
a is a normalizing constant found with a linear fit of the 
background scatter plot constrained to pass through the 
origin. These contours then assign a joint probability to 
the loudest event duple (a^,6^). 

Table II shows sensitivity estimates on GW amplitude 
and GW isotropic emission energy assuming a distance of 
lOkpc to SGR 1900+14, at 90% detection efficiency, for 
the N = 11 flat model and the fluence- weighted model. 
Sensitivity estimates from simulated data (in which there 
could be no GW signal) are produced using the identical 
procedure for producing upper limit estimates in a real 
search. Sensitivity estimates with N = 1 are also shown 
for the sake of comparison; when N = 1 the Stack-a-flare 
pipeline reduces to the individual burst search pipeline 
(Flare pipeline). Superscripts in Table II give a system- 
atic error and uncertainties at 90% confidence. The first 
and second superscripts account for systematic error and 
statistical uncertainty in amplitude and phase of the de- 
tector calibrations, estimated via Monte Carlo simula- 
tions, respectively. The third is a statistical uncertainty 
arising from using a finite number of injected simulations, 
estimated with the bootstrap method using 200 ensem- 
bles [33] . The systematic error and the quadrature sum 
of the statistical uncertainties are added to the final sen- 
sitivity estimates. Methods used to estimate these un- 
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FIG. 10: SGR 1900+14 storm light curve, focusing on fluence (integrated counts) estimates and timescales. The N = 11 flat 
model is illustrated here. In addition to features common to Figure 11, alternating red and cyan highlights indicate integration 
bounds for fluence estimation of the 11 most energetic bursts. Fluence estimation is one of the goals of the light curve 
processing algorithm. The yellow features are 1090 Hz r = 200 ms ringdown waveforms, aligned with the nominal beginning of 
the electromagnetic bursts (indicated as before by green circles). A small constant time offset between this nominal time and 
gravitational wave emission would be immaterial to the search, so long as the time offset were small relative to the on-source 
region half-duration of 2 s. The ±2s on- source regions, which often overlap even in this N = 11 model, are illustrated by the 
green shading. 



certainties are further described in [15]. As mentioned 
in Section IV A, one-sigma burst timing uncertainties as 
measured by fits of burst rising edges are built into the 
Monte Carlo. We present the energy sensitivity estimates 
graphically in Figure 16. 

We also present sensitivity estimates on 7 = 
I -E"em i a measure of the extent to which an energy 
upper limit probes the GW emission efficiency. Estimates 
of 7, being ratios of energies, do not depend on the dis- 
tance to SGR 1900+14. 



V. DISCUSSION AND CONCLUSION 

We have presented the Stack-a-flare method for search- 
ing for GW associated with multiple SGR bursts that ex- 
tends the individual SGR burst search presented in [12, 
15]. We have characterized both the T-Stack and the P- 



Stack versions of the Stack-a-flare pipeline, demonstrat- 
ing sensitivity dependence on stacking number N and un- 
certainty in relative timing between bursts. The P-Stack 
pipeline is robust to timing errors, and we have used it to 
estimate GW amplitude, isotropic emission energy, and 
7 search sensitivities for a mock SGR 1900+14 storm 
multiple SGR search, using simulated data modeled af- 
ter LIGO data from the science run which was ongoing 
at the time of the SGR storm. We considered two GW 
emission models to inform the stacking: flat with N = 11; 
and fluence- weighted with N — 18. 

Two other searches for GWs associated with SGR 
events, besides [12], have been published previously; 
neither claimed detection. The AURIGA collabora- 
tion searched for GW bursts associated with the SGR 
1806—20 giant flare in the band 850-950 Hz with damp- 
ing time 100 ms, setting upper limits on the GW energy 
of ~ 10 49 erg [34]. The LIGO collaboration also published 
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FIG. 11: SGR 1900+14 storm light curve with 1ms bins in the (15-100) keV band. Bottom plot shows a detail. Burst start 
times are estimated by fitting the steeply rising burst edges; EM fluences are estimated by integrating light curve area under 
each burst. A 30-bin running average is shown in addition to the raw light curve. Solid lines are linear fits to rising edges; the 
boundaries of rising edges were found by examining the first derivatives in the neighborhoods of the peak locations. Crosses 
mark burst peaks and intersections of the rising edge fits extrapolated to a linear fit of the noise floor measured in a quiescent 
period of data in the 50 s BAT sequence before the start of the storm. The one-sigma timing uncertainty averaged over all 
measurements is 2.9 ms. X-axis times are relative to 2006-03-29 02:53:09.9 UT at the Swift satellite. 



on the same giant flare, targeting times and frequencies 
of the quasi-periodic oscillations in the flare's x-ray tail 
as well as other frequencies in the detector's band, and 
setting upper limits on GW energy as low as 8 x 10 46 erg 
for quasi-periodic signals lasting tens of seconds [35]. 

Based on the pipeline characterization, we expect to 
gain a factor of y/Tl = 3.3 in energy sensitivity in the 
N = 11 flat model over the N = 1 case. We observe 
a gain of 3.8 averaged over results obtained with the 12 
injection waveforms. Comparing these estimated sensi- 
tivities to the individual burst search results for the SGR 
1900+14 storm published in [12], which analyzed the en- 
tire storm in a single ±20 s on-source region, we observe 
a gain in energy sensitivity of an order of magnitude. 
The best values of 7 in [12], for the 2004 December SGR 
1806-20 giant flare, are in the range 5 x lO 1 ^ x 10 6 de- 
pending on injection waveform. These 7 upper limits are 
about a factor of 10 3 lower than the estimated 7 sen- 
sitivities presented here. However, the Advanced LIGO 
detectors promise an improvement in h ISS by more than a 
factor of 10 over S5, corresponding to an improvement in 



energy sensitivity by more than a factor of 100. A stack- 
ing search similar to this mock search in Advanced LIGO 
data could thus conceivably beat the SGR 1806-20 giant 
flare 7 upper limits. Furthermore, on 2008 August 22, a 
new SGR was discovered [36-38] which may be located at 
a distance of only 1.5 kpc in the direction of the galactic 
anti-center [39, 40]. A stacking analysis on bursts from 
this SGR could therefore gain almost 2 additional orders 
of magnitude in energy and 7 upper limits. 



This method would be suitable for a multiple burst 
search for GW associated with the SGR 1900+14 storm 
using real LIGO data. A real search would be similar 
to the mock search in simulated data presented here, 
although detector collaborations may choose to explore 
different or additional stacking emission models. This 
method would also be suitable for multiple SGR burst 
searches on isolated bursts spanning months or years. 
We expect this method to continue to be useful when 
advanced detectors with greater sensitivity come on line. 
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FIG. 12: Histogram of one-sigma timing uncertainties asso- 
ciated with the rise start time of each burst. These uncer- 
tainties were estimated from linear fits of the rising edges of 
bursts in the light curve assuming errors in the light curve 
data are independent normal with constant variance. The 
smallest uncertainty is 1.6 x 1CP 04 s, the largest is 18 ms, and 
the mean is 2.8 ms. The largest uncertainty is due to a burst 
at ~4.4s which rises out of another large burst, far above the 
noise floor. The other large uncertainties are due to bursts 
with rippling rising edges which may be due to additional 
injections of energy. 
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FIG. 13: Histogram of integrated counts in bursts in the SGR 
1900+14 storm. This quantity is roughly proportional to flu- 
ence. The histogram shows a separation between the 11 most 
electromagnetically energetic bursts and the rest. This sep- 
aration determined the break point for the N = 11 fiat GW 
emission model. 
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FIG. 14: Top plot: cumulative integrated counts. The verti- 
cal line marks the contribution from the 11th most energetic 
peak. 89% of the total integrated counts in the 42 storm 
bursts measured is contained in the most energetic 11 bursts. 
The twelfth most energetic burst contributes an additional 
1% to the total estimated from the 42 bursts. Bottom plot: 
fractional change in the cumulative integrated counts. The 
16th burst contributes less than 1% to the running total. The 
37th burst contributes less than 0.1% to the running total. 
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FIG. 15: Example cumulative histograms showing false alarm 
rates versus analysis event loudness for the background (blue) 
and the stacked on-source region (red). There are three such 
plots for each emission model, one per search band. Since the 
stacked livetime is 4 s here, the loudest on-source event occurs 
once per 4 s, and is plotted at a y- value of 0.25 Hz. We can 
estimate the FAR of this loudest on-source analysis from the 
background. Error bars are Poissonian at 90% confidence. 
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FIG. 16: Stack-a-flare simulated data energy sensitivity es- 
timates at lOkpc, for the 29 March 2006 storm from SGR 
1900+14, for the flat and fluence-weighted models. Uncer- 
tainty estimates have been folded in, as tabulated in Table II. 
Vertical lines indicate boundaries of the three distinct search 
frequency bands. Crosses and circles indicate linearly and 
circularly polarized RDs, respectively. Triangles and squares 
represent 11ms and 100 ms band- and time- limited WNBs, 
respectively. Symbols are placed at the waveform central fre- 
quency. These sensitivity estimates reflect the noise curves of 
the detectors. 
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